CLustering Mouse embryo E8.5

Step 1: Load packages

Script laadt de volgende packages in met de library functie om ze toe te kunnen passen voor analyse:- (dplyr)- (ggplot2) - (pheatmap) - (tidyr)- (RColorBrewer)- (ggrepel) - (Seurat) - (Matrix) - (here)

Step 2: Load Data & create object

Deelvraag: in welk format kan mijn data ingeladen worden?

Data wordt als csv bestand met de “read” functie ingeladen tot object en bewaard in het werkgeheugen. Vervolgens worden de data als object genaamd “counts” ingeladen samen met een matrix tabel (mtx.file). Als laatste wordt dit met behulp van de Seurat packages een object gemaakt als preprocessing voor de data analyse, afkapwaarde is vanaf mininmaal 2000 features (genen).

Step 3: data inspection & quality control

Deelvraag: Welke kwaliteitsstappen heb ik nodig voor een data analyse?

De data moet geinspecteerd worden en gecontroleerd worden op inhoudt. Een manier om dit te doen is met de Seurat package. Deze package heeft ingebouwde functies om de data te controleren op inhoud & data vervuiling (bijvoorbeeld mitochondriaal materiaal).

Dit kan met behulp van het Patroon -mt wordt opgezocht in het object “seurat”. Vervolgens wordt het uitgeplot in een violinplot, net als de functies nFeature_RNA en nCount_RNA. Dit als kwaliteitscheck hoe de data verdeeld is.

## Warning: Default search for "data" layer in "RNA" assay yielded no results;
## utilizing "counts" layer instead.

Violinplot die de verhouding weergeeft van de data aan de hand van de volgende features x-as= datalabel, y-as= aantal : - nFeature_RNA; het aantal genen in de dataset - nCount_RNA; het aantal moleculen in de dataset - percent.mt; het percentage mitochondriaal genexpressie

De data toont een verdeling aan volgens het vioolprincipe, waarbij het begint vanaf 2000 features. Optimaal zou zoveel mogelijk features meenemen beter zijn echter vanwege maximale servercapaciteiten gekozen voor 2000 features om de server te ontlichten en crashes te voorkomen.

Identificatie van inhoud - Weergave het aantal per datalabel

## 
##  E8.5-1 E8.5-10 E8.5-11 E8.5-12 E8.5-13 E8.5-14 E8.5-15 E8.5-16 E8.5-17 E8.5-19 
##     411    1098     778     652      73    1004    1154     865     934     654 
##  E8.5-2 E8.5-20 E8.5-21 E8.5-22 E8.5-23  E8.5-3  E8.5-4  E8.5-5  E8.5-6  E8.5-7 
##     671     938     835     596     798     772     586     696     814     962 
##  E8.5-8  E8.5-9 
##     437     398

Feature Scatter: Inladen van twee nieuwe objecten, waarvan de pearson correlatie geanalyseerd wordt. Dit om aan te tonen of de data in de dataset wel overeenkomt met de verwachtingen, geen data verloren of vervuild is met andere data. De andere genaamd cnt_mt is om aan te tonen dat de set niet vervuild is met mitochondriaal materiaal.

cnt_ftr heeft een correlatie van 0.97 wat aantoont dat het een sterke correlatie heeff, aannemelijk is dat de data zuiver is en niet ontbreekt of vervuild is met andere data.

cnt_mt heeft een correlatie van -0.06 wat aantoont dat het een negatief lage correlatie heeft, we kunnen aannemen dat de dataset weinig vervuild is met mitochondriaal materiaal.

Step 4: Normalize and scale data

Deelvraag: Welke methode van normaliseren kan ik toepassen? Seurat biedt momenteel twee opties aan SCtransform en de standaard LogNormalize, echter vanwege de capaciteit van de server gekozen voor de standaard normalisatie met als vaste parameter 10000.

Data wordt geoptimalizeerd voor betere processing, afkapwaarde zijn vanaf >2000 en maximaal <12000 features, percentaal mitochondriaal is < 5 % Dit ter optimalisatie van de dataset, preprocessing voor de volgende stappen, mindere belasting voor de server en een zo groot mogelijk bereik mee te nemen. Data wordt genormalizeerd met de LogNormalize functie met een factor 10000, hiervoor gekozen omdat dit standaard wordt aangegeven in de Seurat analyse.

## Normalizing layer: counts

Opnieuw gegeneerd een violingplot dit ter aantoning dat de preprocessing goed verlopen is.

Step 5: identifying top genes

Deelvraag: Hoeveel genen neem ik mee voor de PCA analyse?

Het identificeren van genen uit de dataset en een selectie meenemen naar de volgende analyse stap. Selectiemethode is “vst” dat de data op basis van logaritmische variatie standardaardiseert en hiermee een verband schept tussen de data en gen met de hoogste expressie. Gekozen voor een aantal van 2000 hoogste variabelen ter ontlasting van de server en na overleg met opdrachtgever dat dit een standaard keuze is om mee te nemen. Echter betekent dat er wel 42883 geexcludeerd worden, dit vanwege het bereik van de dataset.

## Finding variable features for layer counts

## Warning in scale_x_log10(): log-10 transformation introduced infinite values.

Plot de variabele waarden van de genexpressie. X-as = average expression, y-as is standaard variatie. - De zwarte stippen staan voor de variabelen die niet variabel zijn en dus niet mee zullen genomen naar de volgende processing. Voordeel ontlasting van de server, nadeel minder data mee. Echter vallen deze in de lage expressie radius en willen we vooral kijken naar de hoge expressie. - De rode stippen zijn de variabelen van de 2000 hoogste expressie. Daarvan zijn de 10 met de hoogste expressie weergeven in de plot.

Step 6: scaling genes

Geen deelvraag, wel cruciale stap voor preprocessing PCA analyse. Genen worden op basis van waarden met de geanalyseerde data gescaled. Dit betekent dat de genen die voorkomen in de set aan de data gelinkt worden.

## Centering and scaling data matrix

Step 7: PCA analysis & UMAP

Hoe controleer ik mijn PCA analyse en welke componenten ga ik vergelijken? De data wordt verdeeld in componenten op basis van de grootste variatie middels PCA analyse. De data visualiseren middels UMAP en DIMplot. De controle van de PCA analyse kan middels een SCREE plot en een dimplot om aan te tonen dat data daadwerkelijk in componenten gedeeld is en genen geanalyseerd zijn. (Waarde geven) Gekozen voor 15 componenten na overleg met de opdrachtgever, wegens alleen interesse in een afgeschaalde bereik.

## PC_ 1 
## Positive:  Tafa5-intron, Cadm1-intron, Mllt3-intron, Zbtb16-intron, Dcc-intron 
## Negative:  Dysf-intron, Vim-exon, Igf2-Gm49394-exon, Arhgap31-intron, Hspg2-2210010C04Rik-exon 
## PC_ 2 
## Positive:  Selenop-exon, Sh3tc1-exon, Plxnd1-exon, Cdh5-exon, Stab1-exon 
## Negative:  Ttn-intron, Ttn-exon, Actc1-exon, Slc8a1-exon, Myl4-exon 
## PC_ 3 
## Positive:  Cdh5-exon, Stab1-exon, Fli1-intron, Tie1-exon, Flt1-intron 
## Negative:  Rbp4-exon, Apoa1-exon, Amn-exon, Cubn-exon, Cubn-intron 
## PC_ 4 
## Positive:  Gpc3-intron, Rian-exon, Meg3-Gm27000-exon, Dlk1-exon, Fn1-Olfr1270-Ppbp-exon 
## Negative:  Hba-a1-Hba-a2-exon, Hbb-bh1-exon, Hba-x-exon, Ank1-intron, Ank1-exon 
## PC_ 5 
## Positive:  Mast4-intron, Fhod3-intron, Cdh4-intron, Palld-intron, Nhsl1-intron 
## Negative:  Ahnak-exon, Hba-a1-Hba-a2-exon, Bnc2-intron, Hbb-bh1-exon, Col1a1-exon

Elbowplot: Data opgedeeld in 20 PCA componenten: x-as = componenten aantal, y-as= standaard deviatie (grootse variatie) Print Seurat: Van 5 PCA componentende hoogste positieve waarden en hoogste negatieve waarden weergeven. VizDimloading: x-as= negatieve waarden & PC_component, y-as = naam gen extron/intron

Dit ter controle van de PCA analyse en of het overeenkomt/geen overlappingen zijn.

Visualiseren van de PCA middels UMAP methode, dit om aan te tonen dat data geclusterd is.

Dimplot waarvan x-as= PC_1 en y-as= PC_2. Twee PC componenten uitgeplot tegenover elkaar om overlapping/verschillen aan te tonen.

## Warning: The default method for RunUMAP has changed from calling Python UMAP via reticulate to the R-native UWOT using the cosine metric
## To use Python UMAP via reticulate, set umap.method to 'umap-learn' and metric to 'correlation'
## This message will be shown once per session

Visualisatie van de PCA analyse met UMAP methode, nieuwe berekening laat clusters zien. Bedoeld voor snelle weergave van de verdeling. Legenda van 26 clusters x-as= umap_1 y-as= umap_2; heeft geen definitiefe betekenis, is alleen om een indicatie te geven hoe de clusters zich verhouden. Geen definitieve waarden.